---
title: "Mapping Organizational Culture Schemas Based On Correlational Class Analysis: A Tutorial"
author: 
date: 
output:
  html_document:
    df_print: paged
  code_folding: hide
  pdf_document: default
  word_document: default
---

```{r}
knitr::opts_chunk$set(echo = TRUE)
```


Purpose: In this script file, we detail all the codes needed to run the correlational class analysis (CCA), for more information see readme file. We generate the correlational classes from a set of previously defined variables, illustrating each of the classes with graphs and heatmaps. In next step, we generate clusters with the same variables used in the CCA, with the intention of demonstrating, then that the classes present different results than other grouping techniques. Finally, we illustrate one application of CCA, dividing the sample according to the classes, running regression models for each one, comparing them. 


Loading the necessary packages.

```{r, warning= FALSE, message=FALSE, message=FALSE}
library(tidyverse) # set of packages "tidy"
library(haven) # package to import data
library(janitor) # package to summarize data
library(readr) # to load data
library(dplyr) # to move organized data
library(corclass) # to partition of the data into separate modules
library(ggplot2) # to create and visualize graphics
library(igraph) # to run network graphics and analysis
library(d3heatmap) # to run interactive heatmaps
library(factoextra) # assists in the construction of graphics
library(irr) # to run reliability tests on interrelated variables
library(gmodels) # package with tools for adjusting models
library(psych) # package with quantitative tools
library(jtools) # package for data analysis and presentation
library(stargazer) # package for regression formatting
library(huxtable) # package to assist in formatting tables
library(gap) # we will use to apply the Chow test
```

Loading the database.

```{r, warning= FALSE, message=FALSE}
Data <- read_sav("Organizational_Culture_CCA.sav")
```

Selecting the variables to run the Correlational Class Analysis (CCA). Forming a new object.

```{r, warning= FALSE, message=FALSE}
Data_Class <- select(Data, c(PCoop1,
                                 PCoop2,
                                 PCoop3,
                                 PCoop4,
                                 PCoop6,
                                 PCoop7,
                                 PCoop8,
                                 PCoop9,
                                 PCoop10,
                                 PCoop11,
                                 PCoop12,
                                 PCoop13,
                                 PCoop14,
                                 PCoop15,
                                 PMerit16,
                                 PCoop17,
                                 PMerit18,
                                 PCoop19,
                                 PCoop20,
                                 PMerit22,
                                 PMerit23,
                                 Rig24,
                                 Rig25,
                                 Rig26,
                                 Rig27,
                                 Rig28,
                                 Rig29,
                                 Rig30,
                                 Rig31,
                                 Rig32,
                                 Rig33,
                                 Rig34,
                                 PMerit35,
                                 Rig36,
                                 PMerit40,
                                 PMerit41,
                                 Rig42,
                                 Value45,
                                 Value46,
                                 Value47,
                                 Value48,
                                 Value49,
                                 Value50,
                                 Value52,
                                 Value53,
                                 Value54,
                                 Value55))
```

Running the CCA with 0.05 of significance.

```{r, warning= FALSE, message=FALSE}
cca05 <- cca(Data_Class, filter.significance = TRUE, filter.value = 0.05,
             zero.action = c("drop"), verbose = TRUE)
```

Defining the classes.

```{r, warning= FALSE, message=FALSE}
ccacases05 <- cca05$membership
```

Adding the class variable (CCA) to the database.

```{r, warning= FALSE, message=FALSE}
Data <- Data %>% 
  mutate(CCA = (ccacases05))
```

Plotting the class charts.

```{r, warning= FALSE, message=FALSE}
res1 <- cca05
par(family = 'serif')
plot(res1, 1, cutoff = 0.05, layout.fruchterman.reingold, drop.neg.ties.for.layout = TRUE, main = 'Reactive Class', bw = T)
plot(res1, 2, cutoff = 0.05, layout.fruchterman.reingold, drop.neg.ties.for.layout = TRUE, main = 'Resilient Class', bw = T)
```

Analyzing the correlation of the data that compose the classes. 

Correlation matrix.
```{r, warning= FALSE, message=FALSE}
MatrixCor <- cor(Data_Class, use = "all.obs", method = "pearson")
```

Analyzing the correlation of each class separately.
First it is necessary to make a split in the base, separating the data of class 1 and 2.

```{r, warning= FALSE, message=FALSE}
Data_Class <-  Data_Class %>% 
  mutate(Class = (ccacases05))
Data_Class1 <- subset(Data_Class, Class=="1", -  c(Class))
Data_Class2 <- subset(Data_Class, Class=="2", -  c(Class))
MatrixCor1 <- cor(Data_Class1, use = "all.obs", method = "pearson")
MatrixCor2 <- cor(Data_Class2, use = "all.obs", method = "pearson")
```

Testing the equality of correlation matrices between classes and between the total sample by Jennrich's test (1970).

```{r}
cortest.jennrich(MatrixCor1,MatrixCor2, n1= 105, n2= 102)

cortest.jennrich(MatrixCor1,MatrixCor, n1= 105, n2= 207)
              
cortest.jennrich(MatrixCor2,MatrixCor, n1= 102, n2= 207)
```

Interactive heatmap of the sample correlation matrix, class 1 and class 2. The strongly negative correlations are red and the strongly positive correlations are blue.

```{r, warning= FALSE, message=FALSE}
d3heatmap(MatrixCor1, colors = "RdYlBu", k_row = 2, k_col = 4, xaxis_font_size = 10, yaxis_font_size = 10)
d3heatmap(MatrixCor2, colors = "RdYlBu", k_row = 2, k_col = 4, xaxis_font_size = 10, yaxis_font_size = 10)
```

In order to demonstrate the difference between the CCA and Cluster technique, below we performed a cluster analysis on the same basis and cross-checked with a chi-square test and Kappa to point out the difference between the two.  

Running the Cluster technique with Euclidean distance using the Ward method.

```{r, warning= FALSE, message=FALSE}
Data_Class <-  subset(Data_Class, select =  -  c(Class))
dd <- dist(scale(Data_Class), method = "euclidean")
hc <- hclust(dd, method = "ward.D2")
```

Scree Plot

```{r, warning= FALSE, message=FALSE}
df <- scale(Data_Class)
fviz_nbclust(df, kmeans, method = "wss")  +
geom_vline(xintercept = 2, linetype = 2)
```

Introducing the Dendogram chart.

```{r, warning= FALSE, message=FALSE}
fviz_dend(hc, k = 2,
cex = 0.5,
main = "",
          xlab = "Objects", ylab = "Distance", sub = "",
k_colors = c("#2E9FDF", "#FC4E07"),
color_labels_by_k = TRUE, 
rect = TRUE,
rect_border = c("#2E9FDF", "#FC4E07"),
rect_fill = TRUE,
ggtheme = theme(text=element_text(family="serif", face="bold", size=10))
)
```

Grouping cases into 2 clusters using the K-means method.

```{r, warning= FALSE, message=FALSE}
set.seed(123)
km.res <- kmeans(df, 2, nstart = 25)
```

Defining the cases of each cluster.

```{r, warning= FALSE, message=FALSE}
cluster <- km.res$cluster
```

Adding the cluster variable to the database.

```{r, warning= FALSE, message=FALSE}
Data <- Data %>% 
  mutate(Cluster = (cluster))
```

Cluster plot

```{r, warning= FALSE, message=FALSE}
fviz_cluster(km.res, data = Data_Class,
             axes = c(1, 2),
             main = " ",
palette = c("#E7B800", "#FC4E07"),
ellipse.type = "norm",
ellipse.alpha = 0.5,
star.plot = TRUE, 
repel = TRUE, 
ggtheme = theme_minimal()
)
```

Presenting a cross table between clusters and classes, followed by the chi-square test that indicates an association between both variables.

```{r, warning= FALSE, message=FALSE}
Crosstab <- with(Data, table(CCA, Cluster))
round(prop.table(Crosstab, margin = 2), 2)
plot(Crosstab)
chisq.test(Crosstab)
```

Selecting the variables representing the cluster and classes to apply the Kappa test that presented a moderate result for the variables classes and cluster.

```{r, warning= FALSE, message=FALSE}
Kapa <- subset(Data, select = c(CCA, Cluster))
kappa2(Kapa, weight = "equal")
```

Selecting variables for Regression Analysis.

```{r, warning= FALSE, message=FALSE}
Data.R <- subset(Data, select = c(Employee_satisfaction_and_appreciation,
                                    Cooperative_professionalism,
                                    Hierarchic_rigidity,
                                    Meritocratic_professionalism,
                                    Company_Age_Log, ISO_14000,
                                    Number_Employees_Log,
                                    CCA)
                                    )
```

Diagnosis of variables, analyze: normality, linearity, multicollinearity, extreme values (outliers) and homoscedasticity.

```{r, warning= FALSE, message=FALSE}
pairs.panels(Data.R, col="red")
```

Split the database according to classes.

```{r, warning= FALSE, message=FALSE}
Data.R_Class1 <- subset(Data.R, CCA == "1")
Data.R_Class2 <- subset(Data.R, CCA == "2")
```

```{r}
fit1 <- lm(Employee_satisfaction_and_appreciation ~ 
             Cooperative_professionalism + 
             Hierarchic_rigidity + 
             Meritocratic_professionalism + 
             Company_Age_Log + 
             Number_Employees_Log + 
             ISO_14000, 
           data = Data.R_Class1)
summary(fit1)
```

```{r}
fit2 <- lm(Employee_satisfaction_and_appreciation ~
             Cooperative_professionalism + 
             Hierarchic_rigidity + 
             Meritocratic_professionalism + 
             Company_Age_Log +  
             Number_Employees_Log + 
             ISO_14000, 
           data = Data.R_Class2)
summary(fit2)
```


```{r}
Data.R_Class1 <- Data.R_Class1[-c(5, 15, 16, 20, 29, 48, 52, 55, 58, 59, 61, 77, 80, 85, 90, 91, 92, 95, 102),]
Data.R_Class2 <- Data.R_Class2[-c(7, 18, 20, 61, 65, 69, 70, 77),]
```

Running multiple regressions for the variables: Employee_satisfaction_and_appreciation, Cooperative_professionalism, Hierarchic_rigidity, Meritocratic_professionalism, Company_Age_Log, Number_Employees_Log, ISO_14000. Class 1 only.

```{r, warning= FALSE, message=FALSE}
fit1 <- lm(Employee_satisfaction_and_appreciation ~ 
             Cooperative_professionalism + 
             Hierarchic_rigidity + 
             Meritocratic_professionalism + 
             Company_Age_Log + 
             Number_Employees_Log + 
             ISO_14000, 
           data = Data.R_Class1)
summary(fit1)
plot(fit1)
```

Running multiple regressions for the variables: Employee_satisfaction_and_appreciation, Cooperative_professionalism, Hierarchic_rigidity, Meritocratic_professionalism, Company_Age_Log, Number_Employees_Log, ISO_14000. Class 2 only.

```{r, warning= FALSE, message=FALSE}
fit2 <- lm(Employee_satisfaction_and_appreciation ~
             Cooperative_professionalism + 
             Hierarchic_rigidity + 
             Meritocratic_professionalism + 
             Company_Age_Log +  
             Number_Employees_Log + 
             ISO_14000, 
           data = Data.R_Class2)
summary(fit2)
plot(fit2)
```

Generating and saving in pdf file, the table with the regression templates found for each class, using the function export_summs( ). You can also use the stargazer function of the package with the same name. To be able to save the file in pdf or doc, you must have the MikTex package installed. If you don't, you can view the table in the console by removing the arguments: "to.file" and "file.name".

```{r, warning= FALSE, message=FALSE}
plot <- plot_summs(fit1, fit2, scale = TRUE, model.names = c("Reactive", "Resilient"))
ractheme=theme_bw()+
  theme(panel.border = element_blank(),
        panel.grid.minor=element_blank(),
        axis.title.y=element_blank(),
        axis.line=element_line(),
        text=element_text(family='serif'),
        axis.text=element_text(size=10),
        axis.title=element_text(size=10),
        legend.text = element_text(size = 10))
plot + ractheme
export_summs(fit1, fit2, error_format = "[{conf.low}, {conf.high}]", to.file = "pdf", file.name = "Results.pdf")
```

```{r, warning= FALSE, message=FALSE, results="hide"}
stargazer(fit1, fit2, type = "html", out = "Results.doc")
```

Test Chow. Analyzing whether the two regressions (Class Reactive and Class Resilient) are equal.

```{r}
Data.R_Class1 <- drop_na(Data.R_Class1)
Data.R_Class2 <- drop_na(Data.R_Class2)
y1 <- subset(Data.R_Class1, select = c(Employee_satisfaction_and_appreciation))
x1 <- subset(Data.R_Class1, select = c(Cooperative_professionalism,
                                         Hierarchic_rigidity,
                                         Meritocratic_professionalism,
                                         Company_Age_Log, Number_Employees_Log,
                                         ISO_14000)
                                        )
y2 <- subset(Data.R_Class2, select = c(Employee_satisfaction_and_appreciation))
x2 <- subset(Data.R_Class2, select = c(Cooperative_professionalism,
                                         Hierarchic_rigidity,
                                         Meritocratic_professionalism,
                                         Company_Age_Log, 
                                         Number_Employees_Log, 
                                         ISO_14000)
                                        )
chow.test.r <- chow.test(y1,x1,y2,x2)
print(chow.test.r)
```



## R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see <http://rmarkdown.rstudio.com>.

When you click the **Knit** button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:


## Including Plots

You can also embed plots, for example:

```{r}
tmp <- installed.packages()
installedpkgs <- as.vector(tmp[is.na(tmp[,"Priority"]), 1])
save(installedpkgs, file="installed_old.rda")
```

Note that the `echo = FALSE` parameter was added to the code chunk to prevent printing of the R code that generated the plot.









